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Abstract 

A thin liquid film covered with an insoluble surfactant in the vicin- 
ity of a first-order phase transition is discussed. Within the lubrication 
approximation we derive two coupled equations to describe the height 
profile of the film and the surfactant density. Thermodynamics of the 
surfactant is incorporated via a Cahn-Hilliard type free-energy functional 
which can be chosen to describe a transition between two stable phases of 
different surfactant density. Within this model, a linear stability analysis 
of stationary homogeneous solutions is performed, and drop formation in 
a film covered with surfactant in the lower density phase is investigated 
numerically in one and two spatial dimensions. 

1 Introduction 

The stability and dynamics of thin liquid films have been of considerable interest 
to both experimental and theoretical research [1-4] . When the thickness of a flat 
liquid film is in the range of ~ 100 ran, it becomes sensitive to interaction with 
its substrate. In a certain range of film thickness determined by the exact form 
of the interaction potential, this will render the film unstable with respect to 
small perturbations and a pattern formation process sets in. Depending on its 
initial height, the film breaks up into droplets, labyrinth-like patterns or arrays 
of holes [4,5]. This process is known as spinodal dewetting. 

Brought onto the surface of a liquid film, an insoluble surfactant, for ex- 
ample an organic molecule with a hydrophilic head group and a hydrophobic 
tail group, alters the surface tension and thereby influences the breakup pro- 
cess. In addition, gradients of surfactant density lead to so-called Marangoni 
convection on the surface resulting in new instabilities like surfactant-induced 
fingering [6, 7]. 

Most surfactants exhibit complicated thermodynamics with several phase 
transitions [8,9]. These affect thin film hydrodynamics via an equation of state, 
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Figure 1: Schematic of a surfactant-laden thin film of water on a substrate. 
Height profile h(x,y,t) indicates the film thickness at location (x, y) and time 
t, whereas 7(21, y, t) describes the surfactant density at the surface above (x, y). 



relating surface tension to surfactant density [10]. In studies of the dynamics 
of surfactant covered thin films, surfactant thermodynamics has so far not been 
paid much attention to, because the hydrodynamics is dominated by Marangoni 
convection rather than by effects of lateral pressure and diffusion [11]. 

However, there is a special focus of experimental research on pattern for- 
mation under conditions close to the so-called main transition in monolay- 
ers of lipids like pulmonary surfactant Dipalmitoylphosphatidylcholinc (DPPC) 
[12-14]. In vicinity of this first-order phase transition parts of the surfactant 
in the liquid-expanded (LE) phase and in the liquid-condensed (LC) phase co- 
exist. In these experiments, a substrate is coated with a lipid monolayer via 
Langmuir-Blodgett transfer, i.e. it is withdrawn from a trough filled with water 
on which a lipid monolayer has been prepared. The observed patterns consist 
of ordered arrays of LE and LC domains, including regular stripes and rectan- 
gles, the formation of which is usually attributed to oscillations of the meniscus 
between the water in the trough and the substrate [14]. A full understanding 
of these phenomena cannot be achieved without understanding the dynamics of 
the film and the surfactant near the main transition. The aim of this letter is 
to outline theoretical description of the evolution of a thin film covered with a 
surfactant undergoing a phase transition. 

The model we are going to present is derived within the lubrication ap- 
proximation [1]. We follow the usual approach [7,15-18] to describe the time 
evolution of the surfactant covered thin film by two coupled partial differential 
equations, describing the height profile of the underlying liquid film and the 
surfactant density. The surfactant phase transition is incorporated by choice of 
a suitable free-energy functional which determines the lateral pressure as well as 
the diffusive flux. We are going to perform a linear stability analysis of station- 
ary homogeneous solutions of the derived equations and investigate the effect of 
the surfactant on drop formation by numerical simulations. 
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2 Lubrication approximation 



We consider a thin liquid film covered with an insoluble surfactant on a moving 
solid substrate (see fig. [1]) . The velocity field of the liquid film can be obtained 
within the lubrication approximation [1] . By this procedure the initially three- 
dimensional flow problem is reduced to an effectively two-dimensional one. The 
liquid film is then described by a height profile h(x,y,t), which indicates the 
local film thickness, and the two-dimensional flow field at the surface u(a;, y, t). 
The surfactant density at the surface above the point (x, y) is described by the 
function 7(2;, The continuity equation of an insoluble surfactant has been 
the subject of considerable discussion [19-21]. Surface geometry is of negligi- 
ble influence in the lubrication regime, leaving us with the nondimensionalized 
conservation law 

r T = -v- [ru + il , (1) 

where I is the diffusive flux of the surfactant. Here we have scaled all quantities 
by characteristic values: 

Y _ ex v _ e V rp euo* 

ho tio ho 

U = ~, T = ^L, I=-L, (2 ) 

"o 7o 7o Mo 

and V := (dx,dy) denotes the nabla operator in non-dimensional coordinates 
X, Y. The dimensionless parameter e = ho/lo defines the ratio of characteristic 
height and length scales of the problem. Neglecting surface forces, the flow field 
at the surface U, subject to a no-slip condition at the moving substrate, is given 

by 

U = -— VP + HeC&- 1 V— - — . (3) 
2 cr u Q 

Here H — h/ho denotes the nondimensionalized film height, Vo stands for the 
substrate velocity, 00 describes the surface tension in absence of any surfactant, 
and Ca = /iiio/co is the capillary number with dynamic viscosity fi. Moreover, 
P is a generalized pressure given by 

P = -e 3 C&~ 1 — V 2 H + U(H). (4) 

(TO 

One can see that P contains, besides the Laplace pressure term ~ V 2 £f, the 
disjoining pressure n(iJ) due to interaction of substrate and liquid. In the 
literature, different expressions for the disjoining pressure have been considered 
(see ref. [1] and references therein for a discussion of possible U(H)). Here, we 
will use the expression 

U(H) = ^L_^L (5) 

with the positive Hamaker constants A n ^ m , where m > n, for example n = 
3, to = 9 for Lennard- Jones potentials. The repulsive short-range interaction 
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prevents a complete dry-off and the substrate is always covered with a thin 
precursor film. Gravity could be accounted for by adding a term GH to P, but 
is not considered here since it plays a minor role in the capillary regime. 

For the height profile of the liquid film, we obtain the standard evolution 
equation [1] 



-^VP+^eCa-V^-^ 
3 2 ao uq 



(6) 



It should be noted, that surfactant thermodynamics affects the system in two 
ways. First of all, the diffusive flux I is determined by the chemical potential of 
the surfactant. Second, the presence of a surfactant alters the surface tension 
a of the liquid film, making it dependent on 7 in a way determined by ehe 
equation of state of the surfactant. These two points are discussed in further 
details in the following section. 



3 Surfactant thermodynamics 

The lateral pressure pi a t of a surfactant is defined by [22] 

= ^0 -Plat (7) ■ (7) 

Experimentally, put (7) is usually obtained from surface tensions measurements 
using a film balance, where the available area per surfactant molecule is adjusted 
with the help of a movable barrier [8,9,22]. The resulting isotherms of material 
exhibiting the first-order LE/LC transition display a behaviour reminiscent of 
a three-dimensional van der Waals gas. It is therefore reasonable to model 
the surfactant thermodynamics close to the main transition by a free energy 
suitable for a two-dimensional analogue of a van der Waals gas. Since the 
surfactant density varies along the surface, we apply a Cahn-Hilliard-type free- 
energy functional [23], allowing for a non-uniform free-energy density: 

m = J d D - 1 x{|(V 7 ) 2 + / hom (7)} , (8) 

where k > is constant and /hom( 7 ) denotes the free-energy density of a homo- 
geneous system with surfactant density 7. Assuming the system to be in local 
thermodynamic equilibrium the corresponding lateral pressure is given by [24] 

p la t(7) = -/(7) + 7M (ch ° m) (7), (9) 
where the chemical potential ^( chom ) is obtained from T by functional derivation: 

(chem) = & = _ kV 2 + ^hom _ (1Q) 

So far, there have been no limitations on the choice of /h om - In the spirit of 
Landau's theory of first-order phase transitions [25] we will now restrict ourselves 
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to free-energy densities that can be approximated sufficiently well by a fourth- 
order polynomial around the critical density "f cr of the main transition. Defining 
7 = 7 — 7 cr we obtain 



/bom (7) = /o + /17 + hi + hi + hi 



:4 



(11) 



Realistic values for the parameters fa can be estimated by fitting eq. to 
experimentally obtained isotherms. Notice, that pi a t has to be matched to the 
measured pressure within the coexistence region with the help of a Maxwell 
construction. 

The diffusive current I in eq. (p} is proportional to the gradient of ^,( chcm ) 
with proportionality constant a [26]. Hence, in nondimcnsionalized form, the 
lateral pressure and the diffusive current can be written as 



lat 



e 2 K 

^norr 



\ (vr) 2 



(r cr + f) 



r cr + r) v-r 

<9-Fhom 

df 



(12) 



e A^VA< (chom) 
cro 



eA 



dT 2 



Here, the dimensionless numbers A and K are defined by 



A = 



holouo ' 



K = 



^7o 

/IqCTO 



(13) 



(14) 



and the nondimcnsionalized free-energy density of the homogeneous system is 
given by 



Fho 



where F„ 



CO 



(15) 



Inserting eqs. (fT2|) and (fl~3|) into the evolution equations |T]) and ((6|), we obtain 
the complete set of governing equations: 



H T = - V- 



ff 2 



pea- 1 (1 - P lat (r)) V 2 tf - H(H)} 



-e^-Ca-^flatfr) - -ff— 
2 uo 



(16) 



r T = - v- r 



^V{e 3 Ca- 1 (l-P lat (r))V 2 J ff 



-IL{H)} - effCa-^flatCT) - — 

"o 



e Wr-^vr 
<9r 2 



(17) 
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In the next section, we will investigate the linear stability of stationary homo- 
geneous solutions of these equations. 



4 Linear stability analysis 



In the following we concentrate on the case of a substrate at rest, Vo = 0. For the 
sake of simplicity, we first consider only one-dimensional fields H(X, T), r(X, T). 
Homogeneous film heights and surfactant densities H = H = const., T = T = 
const, are always stationary solutions of the equations. Expanding cqs. (|16p 
and ([H]) like H{X,T) = H + rj{X,T) and T{X,T) = f + C(X,T) yields the 
linearised set of equations 



dn 



A 



(18) 



with the linear operator 



A = 



3 
.H 2 



r- 



dli 
~dH 
dU 
~dH 



0% 



e 3 Ca- 



H 



e 3 Ca- 1 tdj c 



H 



)'( 



d 2 F hc 



dT 2 



d% - e 2 Kd 



Ar 



T Z HG& 



A 



d -Fhoni 



dr 2 



d 2 x 



Ar 



where S = a(^ Q T) / cr and Ar = T — T cr . Using the ansatz r\ ~ exp(AT - 



(19) 
■ikX), 
-A(k) 



^ ~ exp(AT + ifcX), trace r and determinant A of the resulting matrix A 
can be obtained as polynomials of wavenumber k. The growth rates A, which 
are the eigenvalues of A(k), can be calculated by use of 



4A 



(20) 



The system is linearly stable if and only if the conditions r(fc) < and A(fc) > 
are simultaneously fulfilled for all wavenumbers k. Assuming the parameters 
e, Ca _1 , A, K to be positive and taking into account that E is a surface tension 
and hence positive as well, the stability condition can be shown to be equivalent 
to: 

on 



OH 

d 2 F ho , 



H 



dT 2 



> 0. 



> 0. 



(21) 



Ar 



In the following analysis, we will use Pi := dU(H) /dH and P 2 := <9 2 F hom (Ar) /dT 2 
as control parameters. Besides the significance of the sign of Pi, which is well 
known from investigations of spinodal dewetting, there is a similar dependence 
on the sign of P2. This is due to the fact, that a homogeneous distribution 
of surfactant in the spinodal region, where P2 < 0, becomes unstable to spin- 
odal decomposition. Although the stability borders are exactly as would be 
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Figure 2: Stability diagram of stationary homogeneous solutions in the Pi-P 2 
plane. If P\ and P 2 are positive, the solution is stable. Otherwise, there is a 
band of unstable modes < k < k c . 



expected from the isolated subsystems H and T, film and surfactant do not 
decouple linearly as will be discussed in the next section. 

An elementary calculation reveals that whenever condition (|21|) is violated, 
there is a band of unstable modes < k < k c , where growth rate Re(A+) is pos- 
itive. If both, Pi and Pi are negative, there will also be a band of wavenumbers 
with positive Re(A_) reaching from k = to a maximal wavenumber smaller 
than k c . 

In principle, the wavenumber fc max corresponding to maximal growth rate 
Re(A+(A; max )) can be calculated from eq. ([!§)) ■ but for the general case the result 
cannot be stated in a concise manner. However, the upper bound of the band 
of unstable modes, k c , can be calculated analytically and the result depends on 
the signs of Pi and Pi. Defining 



we obtain 



( ki for Pi < 0, P 2 > , 

fc c = < k 2 for Pi > 0, P 2 < , (23) 

[ max{/ci,fc 2 } for Pi<0,P 2 <0. 

This means, that in its lower-left quadrant, the Pi-P 2 plane is divided by the 
line ki — k 2 , or equivalently P 2 = (K /eCa -1 E)Pi, into one region where k c = ki 
and another one, where k c = k 2 (see fig-d]). 

Since operator A contains only even powers of V, it is clear, that by writ- 
ing k := |k|, the same results for A(fc), r(fc) and k c are obtained in the two- 
dimensional case, using the ansatz n ~ exp(AT + ik • X), £ ~ exp(AT + ik • X). 



7 



-0.09 





Plat ( 




I 
i 
i 
i 






0.2 - 




i 
i 
i 

i i 
I i 
1 i 






0.0 - 




1 

1 





0.5 r «' 1.3 r 0.6 r « 1.2 r 1 

Figure 3: The free-energy density (left) used in the numerical example and the 
corresponding pressure-area diagram (right). The solid horizontal line results 
from Maxwell construction. 



5 Numerical analysis 

We have numerically simulated the nonlinear set of equations (fT6)) and (fiTf on 
periodic domains in one and two dimensions using a pseudospectral method of 
lines code [27]. Time integration was performed by an embedded 4(5) Runge- 
Kutta scheme, using the Cash-Karp parameter set [28], while the r.h.s. of the 
evolution equations was calculated using 256 Fourier modes in ID or 64 x 64 
modes in 2D. 

In the simulations we used a simple symmetric double well potential -FLom 
employing the parameters r cr = 1, F x = F 3 = 0, F a = -0.1, F 2 = -0.24, F A = 3.85 
Figure [31 shows the free-energy density Fh om , as is obtained for our choice of pa- 
rameters, as well as the resulting pressure-area diagram calculated according to 
eq. (jT3J) for homogeneous values of T. The two minima of i^om correspond to 
two thermodynamically stable phases of different surfactant density. Within this 
simple model, the phases of higher and lower density can be identified with the 
liquid-condensed (LC) phase and the liquid-expanded (LE) phase, respectively. 

For the disjoining pressure ([5]) we use the parameters n = 3, m = 9 with 
A3 = 3 and A$ = 1 and set the remaining dimensionless numbers to Ca _1 = 
1,A= 0.05, A' = 0.05, v = 0. 

Our goal is to investigate how the surfactant affects the formation of droplets. 
Therefore, as initial conditions, small random perturbations of the stationary 
homogeneous solution H = 1.2 and T — 0.82345 are used. Since 911(1. 2)/dH < 
this corresponds to a flat film, which would be unstable to droplet formation 
even in absence of any surfactant. The value of T is chosen as the position of the 
lower density minimum of -Fhom- Thus, we are simulating a thin film uniformly 
covered with surfactant in its lower density phase. 

The parameter values specified above correspond to the case Pi < 0, P2 > 
of the linear stability analysis. The wavenumber dependent growth rate 
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Figure 4: Eigenvalue spectra of stability matrix A, obtained for the parameter 
set used in our numerical example. 





Figure 5: Space-time plots of the ID liquid film H (X, T) (a) and the surfactant 
density T(X,T) (b), displaying regions of positive and negative curvature in 
black and white, respectively. 
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Re(A+(/c)) is displayed in fig. 0J For comparison we also show Re(A+(fc)) for 
two different sets of H,T. The first one where Pi < 0,P2 < corresponds to 
the unstable fixed point T = T cl = 1 of Fhom and the same H = 1.2, whereas 
the other describes the stable case Pi > 0, P2 > where again T = 0.82345 and 
H = 0.9. Now we are in position to determine numerically the maximal growth 
rates A+(A: max ) and the corresponding eigenvectors for the three parameter sets 
mentioned above. The calculations show that neither component of the eigen- 
vectors is dominant. This indicates the coupling of the fields H and T even 
within the scope of linear theory. 

Now we concentrate again on the initial condition H = 1.2, T = 0.82345. 
Time evolution of the nonlinear system in the one-dimensional case can be de- 
scribed as follows. In the beginning, spinodal dewetting and surfactant spinodal 
decomposition lead to rapid formation of liquid droplets and surfactant domains, 
the latter consisting of surfactant in a higher-density (LC) and a lower-density 
(LE) phase. Then a coarsening process sets in and drops of liquid coalesce 
while surfactant domains merge into larger ones. To visualize the coarsening 
of the liquid film and the surfactant density, we display the regions of negative 
curvature, which naturally indicate the location of drops and domains of high 
surfactant density, in a space-time diagram (see fig. [5]). Figure [5] shows a snap- 
shot of a later stage of the ID simulation, where the few remaining high-density 
domains are clearly located on the largest liquid drops. Evidently, there is a 
strong correlation of T and H . In the 2D case, morphology is similar. As can 
be seen in fig. [3 high-density domains are again located on the largest drops of 
liquid. To emphasise this correlation the contour line r = 1 is drawn in the plot 
of the field H (fig. [7] (a)) and vice versa. Like in the ID simulation, drops, that 
are covered with surfactant in the low-density phase become smaller and smaller 
as T increases, while drops covered with the high-density surfactant phase grow 
stronger during the coarsening. Obviously the surfactant has a sustaining effect 
on the drops, since the system energetically favours surfaces with lower surface 
tension. 

6 Conclusion and outlook 

We have modelled the dynamics of a thin liquid film covered with an insolu- 
ble surfactant in the vicinity of a phase transition. For that purpose we have 
incorporated a suitable free-energy functional for the surfactant into the two 
governing equations, which were derived within the lubrication approximation. 
Linear stability analysis revealed the interplay of surfactant spinodal decompo- 
sition and spinodal dewetting of the liquid film. Although result (|23|) for k c 
seems to imply that the liquid film and the surfactant are linearly decoupled, 
their time evolution is connected from onset. One- and two-dimensional simula- 
tions were presented, showing the decomposition of the surfactant into domains 
of material in thermodynamically stable phases. Droplets and domains show a 
strong spatial correlation resulting from the sustaining effect of the surfactant 
on droplets. The proposed model might serve as a starting point for further re- 
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Figure 6: One-dimensional film profile (solid) and surfactant density (dashed) 
at T = 100. 
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Figure 7: Film height (a) and surfactant density (b), each with contour lines 
(H, T = 1) of the other field, at T = 100. 
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search on thin film dynamics with surfactant phase transitions and the related 
pattern formation. Also, further investigation is needed, to understand the role 
of the phase transition in Langmuir-Blodgett transfer systems. Thus, it will be 
necessary to solve the set of evolution equations subject to suitable boundary 
conditions. 
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